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ABSTRACT 

The transits of a planet on a Keplerian orbit occur at time intervals exactly equal 
to the period of the orbit. If a second planet is introduced the orbit is not Keplerian 
and the transits are no longer exactly periodic. We compute the magnitude of these 
variations in the timing of the transits, St. We investigate analytically several limiting 
cases: (i) interior perturbing planets with much smaller periods; (ii) exterior perturb- 
ing planets on eccentric orbits with much larger periods; (iii) both planets on circular 
orbits with arbitrary period ratio but not in resonance; and (iv) planets on initially 
circular orbits locked in resonance. Using subscript "out" and "in" for the exterior and 
interior planets, fi for planet to star mass ratio and the standard notation for orbital 
elements, our findings in these cases are as follows: (i) Planet-planet perturbations are 
negligible. The main effect is the wobble of the star due to the inner planet, therefore 
St ~ l^iiniain/ Clout) Pout- (h) The exterior planet changes the period of the interior one 
by fJ'outiO'in/'foutYPin- As the distaucc of the exterior planet changes due to its eccen- 
tricity, the inner planet's period changes. Deviations in its transit timing accumulates 
over the period of the outer planet, therefore 5t ^ fJ-outeout{<iout / ciin)^ Pout- (hi) Half 
way between resonances the perturbations are small, of order ^oMtaf„/(ain — ctout)^ Pin 
for the inner planet (switch "out" and "in" for the outer planet). This increases as 
one gets closer to a resonance, (iv) This is perhaps the most interesting case since 
some systems are known to be in resonances and the perturbations are the largest. 
As long as the perturber is more massive than the transiting planet, the timing vari- 
ations would be of order of the period regardless of the perturber mass! For lighter 
perturbers, we show that the timing variations are smaller than the period by the 
perturber to transiting planet mass ratio. An earth mass planet in 2:1 resonance with 
a 3-day period transiting planet (e.g. HD 209458b) would cause timing variations of 
order 3 minutes, which would be accumulated over a year. These are easily detectable 
with current ground-based measurements. 

For the case of both planets on eccentric orbits, we compute numerically the transit 
timing variations for several cases of known multiplanet systems assuming they were 
edge-on. Transit timing measurements may be used to constrain the masses and radii 
of the planetary system and, when combined with radial velocity measurements, to 
break the degeneracy between mass and radius of the host star. 

Key words: planetary systems; eclipses 



1 INTRODUCTION 

The recent discovery of planets orbiting other stars ( "exoplanets" ) has opened a new field of astronomy with the potential 
to address fundamental questions about our own solar system which we can now compare with other planetary systems. The 
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primary mode for discovery of exoplanets has been the measurement of the stellar radial velocities via the D oppler effect. 
Currently the small reflex motion of the star du e to the orbiting planet can only be detected for mpia„et lOm^ fflutler et alJ 
booi iMcArthur et alJbood: ISantos et allbooJ . More recent ly, a large number of planetary transit searches are being carried 
out which are start i ng to yield an handful of giant planets ( Charbonneau et alJl200nl : iKonacki et alJlioO^t IPont et ahlbooj : 



iKonacki et allbooi iBouchv et alJboni lAlonso et aljEooi ) , and many more planned searches should reap a large harvest of 
transiting planets in the near future ijHornel Despite these successes, the discovery of "terrestrial" exoplanets, similar 

in size to the earth, awaits the development of several other techniques such as astrometry, space-based transit searches, 
microlensing, or direct imaging JPerrvmaJ2000l: Icharbonnealboosl : iFord fc Tremain dbood iBorucki et aUbood ICould et alJ 

booJ . 

T he first transiting planetary system, HD 209458b, was discovered with Doppler motions of the primary star JCharb onneau et alJ 
I2OO0I) . Since the mass of the planet is degenerate with orbital inclination, the planetary status of the companion was confirmed 
since the transits imply it is edge-on. HST observations yielded precision measurements of the transit lightcurve (iBrown et alj 
l200ltl . which made this the surest planetary candidate around a main sequenc e star (other than ou r own). The ratio of the 
planetary radius to the stellar radius can be measured with extreme prec ision jMandel fc Aeollboo j). However, the absolute 
radii are uncertain due to a degeneracy between radius and mass of star iSeaeer fc Maflen-Qrnelasll2003ll : an increase in the 
mass and radius of the star can yield an identical lightcurve and period. 

This mass-radius degeneracy may be broken if there is an additional planet in the system. About 10 per cent of the stars 
with known planetary companions have more tha n one planet, while possibly as much as 50 per cent of them show a trend 
in radial velocity indicative of additional planets ijFischer et alJboOiL If one or both of the plane t s is transiting, dynamical 
interactions between t he planets will alter the timing of the transits ([Dobrovolskis fc Boruckilll99fit Icharbonneau et al.ll200fll : 
lMiralda-Escuddl2002l) . A measurement of these timing variations, combined with radial velocity data, can break the mass- 
radius degeneracy. 

Given the dual motivations of detecting terrestrial planets and breaking the mass-radius degeneracy, we derive analytic and 
numerical results for transit timing variations due to the presence of a second planet. We begin our discussion by introducing 
the three-body system in §2. The signal from non-interacting planets is calculated first (§3) and then we compute the effects 
of an eccentric exterior perturbing planet with a large period in §4. A derivation of the general transit timing differences for 
two planets with circular, co-planar orbits is presented in §5. The case of exact mean-motion resonance is analyzed in §6. 
The case of two eccentric planets is considered in §7, along with numerical simulations of several known multi-planet systems 
(these are not transiting). We show how measurements of the dispersion of transit timings can be used to detect a secondary 
planet in the system (§8.1), we compare to other planet-search techniques (§8.2), and we show how to determine the absolute 
size and mass of the objects in the system (§8.3). Finally, we discuss other effects we have ignored that an observer should be 
conscious of (§9). 

Throughout the rest of the paper we characterize the strength of transit timing variations as follows. For a series of transit 
times, tj, we fit the times assuming a constant period, P. We compute the standard deviation, a, of the difference between 
the nominal and actual times. Mathematically, 



1 



^(tj-fo + Pjf 



Lj=o 



1/2 



(1) 



where P and to are chosen to minimize a. If the variations are strictly periodic, then the amplitude of the timing deviation 
is simply y/2 times larger than a. 

During the preparation of this paper a proceedings contribution has appeared by Jean Schneider which considers several 
of the effects discussed here (Schneider 2003); however, we find that Schneider's results are incorrect as he does not consider 
the differential force between the star and the transiting planet. In addition, calculations similar to those presented here are 
being carried out by Matt Holman and Norm Murray. 



2 EQUATIONS OF MOTION 

We are studying the 3-body system in which the three bodies have labels 0, 1, 2 and positions Ri, i = 0, 1, 2 (with an arbitrary 
origin). The exact Newtonian equations of motion are given by 

Multiplying the equations for each particle by its mass and adding together, one finds: 
moRo + miRi 4- 1712^2 = 0. (3) 



© 0000 RAS, MNRAS 000, 000-000 



Transit timing variations 3 




early late on time 

Figure 1. Cartoon showing changes in the timing of transit due to a perturbing planet interior to the orbit of a transiting planet. 



This is simply a statement that the centre of mass of the system has no external forces. Since hght travel time and parallax 
effects are negUgible ( H9.1II . the transit problem is unaffected by the total velocity or position of the centre of mass, so we set 

Rc™ = ^? = 0. (4) 

This reduces the differential equations of motion to two, which we take to be that of the two planets, Ri and R2 (for the two 
planetary masses). We use this system of equations for numerically solving the equations of motion. However, for analytic 

consideration it is more convenient to write the problem in Jacobian coordinates which we discuss next. 

The Jacobian coor dinate system is commonly used in perturbation theory for many bodies (see, e.g. iMurrav DermottI 
11994 lMalhotralll993allbl) . For the three-body problem, the Jacobian coordinates amount to three new coordinates which 
describe (a) the centre of mass of the system; (b) the relative position of inner planet and the star (the "inner binary" ) ; (c) 
the relative position of the outer planet and the barycentre of the inner binary (the "outer binary"). To distinguish from the 
body coordinates, we denote the Jacobian coordinates with a lower case r^. The Jacobian coordinates are 

ro = R-cm ~ 0, 
ri = Ri — Ro, 

r2 = R2 ; [moRo + miRi] . (5) 

mo + mi 

Using — rrii/M ~ m^/mo, where M = X^^^q rrii, the equations of motion may be rewritten in Jacobian coordinates, 

Gmo ri ri-r2i r2i 

ri = --GM^j.2-. p7-GM/i2^, 

1 — ^1 rj jn — r2i|'^ 

Gmo r2i V21 - ri , , 

V2 = 3--GM/ii| -, (6) 

1 - '■21 |r2i - ri|3 

where V21 = /iiri + r2 = R2 — Ro- 
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3 NON-INTERACTING PLANETS: PERTURBATIONS DUE TO INTERIOR PLANET ON A SMALL 
ORBIT 

Throughout the rest of the paper we make the approximations that (a) the orbits of both planets are ahgned in the same 
plane; (b) the system is exactly edge-on, that is, the inclination angle is 90°. We also approximate the planet and star as 
spherical so that the transit is symmetric with a well-defined midpoint. 

If we take the limit as /ii,^2 — * in equation ©, then the orbits of the planets follow Keplerian trajectories with the 
equations of motion 

r2 = -Gmo^. (7) 

'2 

This approximation requires that the periapse of the outer planet be much larger than the apoapse of the inner planet, 
(1 — 62)02 ^ (1 -|- ei)ai where ai, 02 are the semi-major axes of the inner and outer binary and ei, 62 are the eccentricities. 
In this case, the inner binary orbits about its barycentre which in turn orbits about the barycentre of the outer binary but 
there is no perturbation to the relative motion of the inner binary due to gravitational interactions. Timing variations that 
arise are simply due to the reflex motion of the star (as shown in Figure 0. 

The simplest case to consider is that in which both the inner and outer binary are on approximately circular orbits. 
The transit occurs when the outer planet is nearly aligned with the barycentre of the inner binary and its motion during the 
transit is essentially transverse to the line of sight. The inner planet displaces the star from the barycentre of the inner binary 
by an amount 

a;o = -ai^i sin [27r(t - to)/A], (8) 

where the inner binary undergoes a transit at time to and Pi is the orbital period of the inner binary. Thus, the timing 
deviation of the mth transit of the outer planet is 

0, xo P2ai/xi sin [27r(mP2 - to)/A] 

dt2 « ~ , (9j 

V2 — vo Z-Ka2 

where Vi is the velocity of the ith body with respect to the line of sight. Typically «o ^ «2, so we have neglected vo in the 
second expression in the previous equation. 

Computing the standard deviation of timing variations over many orbits gives 

IHi-^ \2\l/2 P2ffllMl .,„s 
Cr2 = {{5t2) ) = . (10) 

Note that if the periods have a ratio P2-P\ of the form q:l for some integer q, then the perturbations disappear since the 
argument of the sine function is the same each orbit. Another observable is the duration of the transit, which scales as 

t2-^^. (11) 

(V2 — VQ ) 

This leads to significant variations only if iiq — W2, or /i^ > 01/02, which requires a very large axis ratio. 

More interesting variations occur if either or both planets are on eccentric orbits. Because both planets are following 
approximately Keplerian orbits, the transit timing variations and duration variations can be computed by solving the Kepler 
problem for each Jacobian coordinate. Since we are assuming that the planets are coplanar and edge-on, 4 coordinates each 
suffice to determine the planetary positions: ei.2, 01.2, vj\^2 (longitude of pericentre), and /i,2 (true anomaly). As in the 
circular case, the change in the transit timing is approximately 5t2 ~ X(i/v2- The position of the star with respect to the 
barycentre of the inner binary is 

a:o = -A*i'"i sin [/i ii7i]. (12) 

If oi <C 02, the outer planet is in nearly the same position at the time of each transit and its velocity perpendicular to the 
line of sight is 

2-7ra2 (1 -f e2Costi72) 

V2 = , , (13) 

P2y/l^2 

where we have used the fact that /2 = ~txj2 at the timing of the transit. Thus, to first order in ai/02 



_ P2Mi'-i sin [/i + roi] ^1 - e| 

27r02(l-|-e2C0Sti72) 

The standard deviation of it2 can be computed analytically as well. Over many transits by the outer planet, the inner 
binary's position populates all of its orbit provided the planets do not have a period ratio that is the ratio of two integers. 
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Consequently, we find the mean transit deviation by averaging over the probabihty that the inner binary is at any position in 
its orbit, p(/i) — (where ni — 27r/Pi), times the transit deviation at that point. This gives 



3 a\ 

-fii — ei smci7i. 

2 V2 



(15) 



Since the star spends more time near apoapse, the mean timing grows as ei. The symmetry of the orbit about zu = and tt 
explains the dependence on sinT^Ji. A similar calculation gives ((Jt^) ^nd the resulting standard deviation is 



(^2 



P2aifj,iy^l — e 



£1 

2 



(1 + 



cos ZUi I 



1/2 



(16) 



2^/27102(1 + 62 COS 1*72) 

This agrees with equation IIOII in the limit ei 0. Averaging again over zui and VJ2, gives 

(<T2)^i.^, 23/27ra2(l-ei)V4 

Note that an eccentric inner orbit reduces (J2 because the inner binary spends more time near apoapse as the eccentricity 
increases thus reducing the variation in position when averaged over time. As ei approaches unity for an orbit viewed along 
the major axis, (T2 reduces to zero since the minor axis approaches zero leaving no variation in the xo position. 



(17) 



4 PERTURBATIONS DUE TO EXTERIOR PLANET ON A LARGE ECCENTRIC ORBIT 

In this section we include planet-planet interactions and compute the timing variations due to the presence of a perturbing 
planet on an eccentric orbit with a semi-major axis much larger than that of a transiting planet on a nearly circular orbit. 
In this limit, resonances are not important and the ratio of the semi-ma jor axes can be used as a small parameter for a 
perturbation expansion. A general formula for this case has been derived bv lBorkovits et all (l2003f) . Here we present a shorter 
derivation which clarifies the primary physical effects for coplanar planets viewed edge-on. 

The equations describing the inner binary can be divided into a Keplerian equation lO and a perturbing force proportional 
to 7712. The perturbing acceleration Sr'i on the inner binary due to the outer planet is given by 

Sri = -GMfi2 ~ ^^].. - GM^j,2^. (18) 



We expand this in a Legendre series and keep terms up to first order in the ratio of the radii of the inner and outer orbit, 

+ 0{ri/r2f. (19) 



J .. Gm; 
ori = =- 



ofi • r2 
ri - 3 5 — r2 



To compute the perturbed orbital period we must find the change in the force on the inner binary due to the outer planet 
averaged over the orbital period of the inner binary. Since the inner binary is nearly circular, the angle of the inner binary is 
given by 9i — f\ + tui ~ ni(t — t\)zu\, where we have approximated ei ~ 0. Differentiating this with respect to time gives 

Oi = hi{t — Ti) + ni — n-iTi. (20) 

Now, we write ni = — 3ni/( 2ai)ai, and e xpress ai, zui, and fi in terms of the radial, tangential, and normal components of 
the force (see section 2.9 of lMurrav fc Dermott,,1999,'l . Plugging these expressions into 9\ gives a cancellation of most terms 
to lowest order in ei, and after setting the normal force to zero leaves the remaining term 

^1 = ^1 fl - 7v^^-^^ ' (21) 
\^ G(mo-fmi) J 

where R is the radial disturbing force per unit mass, R = ((5fi • fi) — ^Gm2ai/r2- 

Thus, the presence of the second planet causes a change in the effective mass of the inner binary by an amount 
— ■^'rn2{ai/r2)'^, which results in a slight increase in the period of the orbit. The increase in period would be constant if 
the second planet were on a circular orbit. However, for an eccentric orbit, the time variation of r2 induces a periodic change 
in the orbital frequency of the inner binary with period equal to P2. 

Now, the time of the (A'^ + l)th transit occurs at 

j.fo + 2TTN 
tea — to = / rf/l^l 
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Figure 2. Contour plot of the logarithm of the dispersion of the normalized timing variations, log {(Tini/x^^), for an inner circular planet 
and an outer eccentric planet (for example, at the -2 contour an orbit lasting 27r years with a perturbing planet of mass 1O~^M0 would 
have a transit timing standard deviation of 10~^ years, or 5 minutes). The dotted line is the approximation given in equation 1271 . 



/o+27rJV 



1 + 



7712 



(22) 



mo + mi 

where /o is the true anomaly of the inner binary at the time of the first transit. Following iBorkovits et all (I2OO3II . we change 
the variable of integration from /i to /2, the true anomaly of the outer planet, 

P2 



dfi = 



^2 



Piai(l-ei)i/2 



d/2 



(23) 



Since we are assuming that the orbit of the outer planet is eccentric, r2 — 02(1 — e2)/(l + 62 cos /2), which gives the transit 
time 

-3/2 



to = NPi 



m2 ( 1 



Pi 



if 2 +62 sin /a) 



(24) 



27r(mo + mi) 

where /2 is the true anomaly of the outer binary at the timing of the (A'^ — 1) transit. The unperturbed /2 includes the mean 
motion, n2{t — T2), which grows linearly with time. To find the deviation of the time of transits from a uniform period, we 
subtract off this mean motion as well as NPi which results in 

3/2 



Sti 



m2 



[f2 
P? 



- n2{t- T2) + 62 sin/2] 



(25) 



27r(mo + mi) P2 

This agrees with the expression of IBorkovits et all (l2003tl in the limit 1 = 1 (i.e. coplanar or b its). R emarkably, the timing 
variations scale as a^'^'''^, which is a much shallower scaling than estimated bv iMiralda-Escudl j20o3), 02^^- 

Numerical calculation of the 3-body problem show that this approximation works extremely well in the limit ri ^ T2 (see 
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Figure |5J. If Pi <C P2 and the period ratio is non- rational, then over a long time the transits of the inner planet sample the 
entire phase of the outer planet. Thus, we can compute the standard deviation of the transit timing variations as in equation 



27r 



ai = {Stiy/^ = ^ dhSMh), (26) 
since {5t\) — 0, where p(/2) ~ 112/ f2- This integral turns out to be intractable analytically, but an expansion in 62 yields 



3/362 

o-i — 



3 2 47 4 413 6 

— e-i go 69 

16 1296 27648 



1/2 

(27) 



^2(1- 02 

which is accurate to better than 2 per cent for all 62. Figure |5| shows a comparison of this approximation with the exact 
numerical results averaged over izj2 (since there is a slight dependence on the value of 1172). This approximation breaks down 
for 02(1 — 62) < 5ai since resonances and higher order terms contribute strongly when the planets have a close approach. It 
also breaks down for 62 < 0.05 since the perturbations to the semi-major axes caused by interactions of the planets contribute 
more strongly than the tidal terms which become weaker with smaller eccentricity. 



5 PERTURBATIONS FOR TWO NON RESONANT PLANETS ON INITIALLY CIRCULAR ORBITS 

In this section we estimate the amplitude of timing variations for two planets on nearly circular orbits. The resonant forcing 
terms are most important in determining the amplitude, even for non-resonant planets. The planets interact most strongly 
at conjunction, so the perturbing planet causes a radial kick to the transiting planet giving it eccentricity. Since the planets 
are not exactly on resonance, the longitude of conjunction will drift with time, causing the kicks to cancel after the longitude 
drifts by ~ TT in the inertial frame. Thus, the total amplitude of the eccentricity grows over a time equal to half of the period of 
circulation of the longitude of conjunction. The closer the planets are to a resonance, the longer the period of circulation and 
thus the larger the eccentricity grows. The change in eccentricity causes a change in the semi-major axis and mean motion. 

For two planets that are on circular orbits near a j:j + 1 resonance, conjunctions occur every Pconj = 27r/(ni — 722) — jP 
(we take the limit of large j and we ignore factors of order unity). We define the fractional distance from resonance, e = 
|1 — (l+j~^)Pi/P2| < 1, where e = indicates exact resonance. Then, because the planets are not exactly on resonance, 
the longitude of conjunction changes with successive conjunctions and the longitude of conjunction returns to its initial value 
over a period Pcyc = Pj~^^~^- The number of conjunctions per cycle is Nc = Pcyc/ Pconj — 3~^^~^ ■ Each conjunction changes 
the eccentricity of the planets by Ae ~ fip^rtia/ ^ci)^ (using the perturbation equations for eccentricity and the impulse 
approximation, where ^ipert is the planet-star mass ratio of the perturbing planet). Over half a cycle the eccentricities grow 

to about NcAe ~ /ipert(l — A/P2)~^(je)^^ — Mperte"^. 

To find the change in the transit timing, we use the orbital elements to compute the variation in the instantaneous orbital 
frequency, 9. To first order in e 

• n(l + e cos f)'^ or,-, / ^ 

~ no + 5n + 2enocos[\-vj\, (28) 



(1 



x3/2 



where no is the unperturbed mean motion. There are two terms which contribute to timing variations: fluctuations in the 
mean-motion and fluctuations due to a non-zero eccentricity. In the flrst case, 5n may be found by applying the Tisserand 
relation to the lighter planet (we now use subscripts "light" and "heavy"), resulting in Snugut / riught — j{Aeiight)^ — 
(where /i is fiheavy)- These changes to the period accumulate over an entire cycle, giving 

StHght ^ f/e-^P. (29) 

By conservation of energy, the fractional change in semi-major axis (or period) of the heavy planet is reduced by a factor of 
miight/mheavy, SO that 

Stheavy — {miight / mheavy) fi^ P- (30) 

The eccentricity dominated term gives a timing variation of 

St ~ ^pert(~^P- (31) 

So, for e > ^^^'^ the perturbed eccentricity dominates, but closer to resonance for j^^^ fj?^^ < e < jj}''^ the perturbed mean 
motion dominates (this range is the same for both the light and heavy planets, except for factors of order unity). For smaller 
values of e, the planets are trapped in mean-motion resonance, which is discussed in the next section. Half way between 
resonances, e ~ so the timing deviation become 

Str^0.7^Ipert{a/{a2-al)fP. (32) 
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P /P -1 
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Figure 3. Transit timing standard deviation for two planets of mass mtrans = 10~^mo and rripert = 10~^mo on initially circular orbits 
in units of the period of the transiting planet. The solid (dotted) line is the numerical calculation for the inner (outer) planet averaged 
over 100 orbits of the outer planet with the planets initially aligned with the observer. The dashed line is equation ilOI . The large dots 
are equations ^381 on resonance and equation 1321 halfway between resonances. 



A more precise derivation in the eccentricity-dominated case using perturbation theory is given in Appendix El 

So far we have discussed the timing variations for planets nearby a first order resonance. For larger period ratios, the 
eccentricity of the inner planet grows to etn ~ fiout{Pin/Pout)'^ , so Sti„ ~ fiout{Pin/ Pout)'^ Pin- For an outer transiting planet 
the motion of the star dominates over the perturbation due to the inner planet for Pout > (2 

Figure 13 shows a numerical calculation of the standard deviation of the transit timing variations. We have used small 
masses to avoid chaotic behavior since resonant overlap occurs for j > /i"^''^ iWisdomlll980l) . Figure 2] zooms in on the 2:1 
resonance. As predicted, the amplitude scales as (equation I31II . and then steepens to (equations 1301 and I29II closer to 
resonance. Since the strength of the perturbation is independent of whether the perturbing planet is interior or exterior, the 
strength of the resonances are similar and the shape of the standard deviation of the transit timing variations is symmetric 
about Pin = Pout- The dashed curve in Figure |3 shows the analytic approximation from equation HUH , which agrees well for 
Ppert < (27r)~^''*Ptra„3. The numerical results match the perturbation calculation, equations lA7t and lASI l. except for near 
resonance where the change in mean-motion dominates (we have not bothered to overplot the perturbation calculation since 
it is indistinguishable from the numerical results). 

There is a dip in fj 2 iiGcir -I out — 2-5Pin which occurs because the amplitude of the timing differences due to the orbit 
of the star about the barycentre (eqn. IIOI I are opposite in sign and comparable in amplitude to the differences due to the 
perturbation of the outer planet by the inner planet (eqn. IA8t . 

The analysis in this section breaks down near each resonance because we have not considered changes to the orbital 
elements of the perturbing planet. In the next section we consider what happens to planets trapped in resonance. 
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6 PERTURBATIONS FOR TWO PLANETS IN MEAN-MOTION RESONANCE 

The analysis in the previous section assumes that the perturbation to the orbit of each planet is small, so that the interaction 
can be calculated using the unperturbed orbits (linear perturbation theory) . This is clearly not the case near a mean-motion 
resonance. We investlRate the case of low, initially zero, eccentricity where we found the standard analyses of this case (e.g. 
iMurrav fc Dermottlll999l) to be incorrect. Here we provide a physically motivated, order of magnitude, derivation of the 
perturbations and the transit timing variations for two planets in a first-order mean-motion resonance. A rigorous derivation 
is left for elsewhere, but we verify our findings with numerical simulations. 

Consider a first order, j'-j+l, resonance where the lighter planet is a test particle. Qualitatively, the physics of low 
eccentricity resonance is as follows: on the nominal resonance, the two planets have successive conjunctions at exactly the 
same longitude in inertial space. The strong interactions that occur at conjunctions build up the eccentricity of the test 
particle and cause a change in semimajor axis and period. The change in period of the test particle causes the longitude of 
conjunction to drift. Once the longitude of conjunction shifts by about tt relative to the original direction, the eccentricity 
begins to decrease making a libration cycle. The libration of the semi-major axes causes the timing of the transits to change. 

This qualitative discussion leads directly to an estimate of the drifts in transit times. Within each libration cycle the 
longitude of conjunction shifts by about half an orbit, mostly due to the period change of the lighter planet. Since conjunctions 
occur only once every j orbits the largest transit time deviation of the lighter planet during the period of libration is P/j 
(in this order of magnitude derivation we ignore factors of order unity, and take the limit of large j so that j ~ j + 1 and 
P2 — Pi). The observationally more interesting case is probably that in which the heavier planet is the transiting one. Then, 
conservation of energy for the orbiting planets implies that the change in periods is inversely proportional to the masses, 
therefore the timing variations are given by {miight/niheavy)P/ j ■ We find an excellent fit to the data for 

Ttlpert ~r Tfitrans 

The calculations shown in Figure verify this analytic scaling with j. 

Calculating the libration period is a little more complicated, but still straightforward. Suppose the period of the test 
particle deviates from the nominal resonance by a small fraction e. Then, consecutive conjunctions drift in longitude by about 
27rj^e. The number of conjunctions, A^c, before a drift of order tt in the longitude of conjunctions accumulates is Nc ~ j~^e~^ . 
We now estimate e indirectly. The test particle gains an eccentricity of order in each conjunction due to the radial force from 
the massive planet (this can be computed from the impulse approximation and the perturbation equation for eccentricity). 
The eccentricity given in Nc conjunctions is then of order Ae ~ iJLe~^ . Using the Tisserand relation, the fractional change in 
semimajor axis associated with this change in eccentricity is j/i^e~^. Since this is also the fractional change in the period we 
have e ~ j^^^ 1^^^^ and a libration period of 

Pub ~ 0.5j-'e-'P ~ Q.^^^^'■'^l'^'''■P. (34) 

We numerically computed the amplitude and period of the transit timing variations at the 2:1 resonance. Figure 01 shows 
a plot of the amplitude of the timing variations versus the mass ratio of the perturbing planet to the transiting planet. As 
predicted, the amplitude is of order the period of the transiting planet when the transiting planet is lighter, and varies as the 
mass ratio when the transiting planet is heavier. The libration period measured from the numerical simulations shows the 
predicted behavior, scaling precisely as pT"^^^ for the more massive planet (with a coefficient of ~ 0.7 for j — \ and 0.5 for 
j > 1 in equation I34II . We have compared the numerical values of the amplitude and period of libration on resonance as a 
function of j. Despite the fact that the above scalings were derived in the large- j' limit, the agreement is better than 10 per 
cent for j ^ 2, and accurate to about 40 per cent for j = 1. 

Figure0]shows the more detailed behavior of the amplitude near the 2:1 resonance. The amplitude is maximum slightly 
below resonance at the location of the cusp. This may be understood as follows: since the simulations are started with 
d = 62 = 0, after conjunction the eccentricity grows and the outer planet moves outwards, while the inner planet moves 
inward. This causes the planets to move closer to resonance, causing a longer time between conjunctions, leading to a larger 
change in eccentricity and semi-major axis. The cusp is the location where the planets reach exact resonance at the turning 
point of libration, at which point 5t is maximum. To the right of the cusp, the libration causes the planets to overshoot 
the resonance, so the change in eccentricity and semi-major axis is somewhat smaller, and hence the amplitude is smaller. 
Figure |1] shows that the width of the resonance scales as ii'^^ (the horizontal axis has been scaled with ji'"^^^ so that the 
curves overlap), so for larger mass planets the resonant variations have a wider range of influence than the non-resonant 
variations discussed in the previous section. The curves in Figure |1] demonstrate that on-resonance the amplitude scales as 
min{l, Upert/ fJ-trans) /j , whilc off-rcsonance the amplitude scales as /ipert. 
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-30 -20 -10 10 20 30 

(P/P,-2)m-^^^ 

Figure 4. Amplitude near ttie 2:1 resonance versus tlie difference in period from exact resonance for two systems: one with mi = 10~^mo 
(top solid) and m2 = 10~^ (lower solid), and tfie other with mi = 10~®mo (top dotted) and m2 = 10~^ (lower dotted). The large dots 
are equation 1331 . The dashed line is equation 1311 . while the dash-dot line is equation I3UI . 



7 NON-ZERO ECCENTRICITIES 

When either eccentricity is large enough, higher order resonances become important. In particular, the resonances that are 
l:m begin to dominate as the ratio of the semi- major axes beco mes large; as the e ccentricity of the outer planet approaches 
unity these resonances become as strong as first order resonances ijPan fc SarJ20o3) . Figure|^shows the results of a numerical 
calculation where the transiting planet HD 209458b, with a mass of approximately 0.67 Jupiter masses, is perturbed by a 
IM0 planet with various eccentricities (we have taken HD 209458b to have a circular orbit). Near the mean- motion resonances 
the signal is large enough that an earth-mass planet would be detectable with current technology. The amplitude increases 
everywhere with eccentricity. This graph can be applied to systems with other masses and periods as the timing variation 
scales as St oc Ptrans'mpert (except for planets trapped in resonance). 

When both planets have non-zero eccentricity, the parameter space becomes quite large: the 4 phase space coordinates for 
each planet, assuming both are edge-on and 2 mass-ratios give 10 free parameters. On resonance, the analysis remains similar to 
the circular case. The libration amplitude will still be of order ~ Pj~^ for the lighter planet and ~ P {^-ught / t^heavy) for the 
heavier planet. However, the period of libration will decrease significantly as the eccentricity increases since Pm ~ e~"'^^^/i~"^^^. 

O n the secular time-scale, the precession of the orbits will lead to a significant variation in the transit timing f Miralda-Escudd 
l2002fl . The period of precession, P^, may be driven by other planets, by general relativistic effects, or by a non-spherical stel- 
lar potential, but lead s to a magnitude of transit timing deviation which just depends on the eccentricity for P^ ^ P. 
iMiralda-Escudd l)2002ll showed that the maximum deviation for e 1 is given by 

eP 

3t^— (35) 

TT 

and the timing variations vary with a period that is equal to the period of precession. For arbitrary eccentricity, the maximum 
deviation is 
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Figure 5. Dispersion, atrans, of timing variations of HD 209458b due to perturbations induced by an exterior earth-mass planet with 
various eccentricities and periods. The color bar has crtrans for a planet of mass and for a transiting period of 3.5 days. Near 
the period ratio of 4:3 the system becomes chaotic. The increase in the number of resonance peaks, particularly those beyond the 2:1 
resonance, are from higher order terms in the expansion of the Hamiltonian which were truncated for the near-circular case. 



St ■ 



p_ 

2^ 



sin y + sin z + 



(36) 



where x — {1 ~ e^)^''*, y = (1 — x)/e, and z — {1 ~ x'^)/e (this is derived from the Keplerian solution with a slowly varying 
■oo). This approaches P/2 as e 1. Typically the eccentricity will vary on the secular time-scale, so these expressions only 
hold as long as the variation in e is much smaller than its mean value. 

Rather than systematically studying the entire parameter space, we now investigate several specific cases of known 
extrasolar planets to demonstrate that detection of this effect should be possible once a transiting multi-planet system is 
found. Most of these systems have non-zero eccentricities and several are in resonance, causing a significant signal. We 
summarize the amplitude of the signals of most known multi-planet systems, if they were seen edge-on, in Table 1 (in some 
cases other planets are present which would cause additional perturbations). 

The extrasolar planetary system Gliese 876 contains tw o Jupiter-mass plan ets on modestly eccentric orbits which are near 
the 2:1 mean-motion resonance, Pi = 30. Id and P2 = 61. Od ijMarcv et al.ll200lh . Due to the small size of the M4 host star, the 
inner planet has a 1.5 per cent probability of transiting for an observer at arbitrary inclination. The orbital motion involves 
both mean-motion resonance as well as a secular resonance in whi ch the planets librate about th eir apsid al alignment. The 
apsidal alignment is in turn prece ssing at a rate of —41° per year jLaughlin et aliboO^I : iNauenberg 2002: .Rivera fc Lissau'eil 
boOltlLaughlin fc ChamberslboOlh . Figure El shows t he predicted timing va riations if this system were seen edge-on and if the 
planets are coplanar using the orbital elements from iLaughlin et al.l l)200'J) . 

The two most prominent periodici ties in Figure El are those ass ociated with the 2:1 libration, with a period of roughly 
600 days (20 orbits of the inner planet, iLaughlin fc Chamber jEool . and the long term precession of the apsidal angle with 
a period of about 3200 days (110 orbits of the inner planet, corresponding to —41° yr~^). Evaluating equation 13511 gives a 
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Table 1. Timing variations for known multi-planet systems 



System 


Prn (d) 


Pout / Pin 




0-2 


55 One e, b 


9 SI 
Z.Ol 


^ 91 


iU.O s 


Z.DO S 


55 Cnc b, c 


14.7 


3.02 


1.61 h 


14.7 ii 


Ups And b, c 


4.62 


52.3 


1.30 s 


1.61 min 


Gliese 876 


30.1 


2.027 


1.87 d 


14.6 ii 


HD 74156 


51.6 


39.2 


4.98 min 


42.4 min 


HD 168443 


58.1 


29.9 


12.9 min 


2.62 ii 


HD 37124 


152 


9.81 


3.43 d 


11.2 d 


HD 82943 


222 


2.00 


34.9 d 


30.7 d 


PSR 1257+12 b, c 


66.5 


1.48 


15.2 min 


22.3 min 


Earth/ Jupiter 


365 


11.9 


1.42 min 


24.1 s 



peak amplitude of 1.4 days for the inner planet and 18 hours for the outer planet which both compare well with the numerical 
results given that the eccentricities are not constant. 

The extrasolar planetary system 55 Cancri contains a set of planets, h and c, near the 3:1 resonance having 15 and 45 
day periods. There is some evidence for another planet, d, i n an extremely long orbit, and recently a fourth low mass planet, 
e, was found with a 2.8 day period (|McArthxir et al]l20n4h . The planets e, 6, and c have transit probabilities of 12, 4, and 
2 per cent, respectively, for an observer at arbitrary inclination. The orbit of planet h is approximately circular while planet 
c is somewhat eccentric jMarcv et al.ll2002l) . Table 1 gives the amplitude of the variations for the planets. We have ignored 
planet e; however, it is at a large enough semi-major axis to produce a ~ 22 second variation due to light-travel time as the 
barycenter of the inner binary orbits the barycenter of the triple system were the inner planets transiting. 

The double planet system Upsilon Andromedae has a semi-major axis ratio of 14 which is not in a mean-motion resonance 
JSutler et alJll999l: iMarcv et aljl20ol . The inner planet has a short period of 4.6 days, and thus a significant probability of 
transiting of about 12 per cent, but has variations which are too small to currently be detected from the ground or space. The 
outer planet has much larger transit timing variations due to its smaller velocity, but a much smaller probability of transiting. 

The planetary sys tem HD 37124 has two planets with a period ratio of ~ 10 and a period of the inner planet of 241 
days jVogt et aljl200(l) . The outer planet is highly eccentric, 62 = 0.69, and so its periapse passage produces a large and 
rapid change in the transit timing of the inner planet. If this system were transiting, the variations would be large enough 
to be detected from the ground. HD 82943 is in a 2:1 resonance giving variations of order the periods of the planets. The 
pulsar planets are near a 3:2 resonance, which would cause large transit timing variations were they seen to transit the pulsar 
progenitor star. Finally, alien civilizations observing transits of the Sun by Jupiter would have to have ~ 10 second accuracy 
to detect the effect of the Earth. 



8 APPLICATIONS 

8.1 Detection of terrestrial planets 

The possibility of detecting terrestrial planets using the transit timing technique clearly depends strongly on (1) the period of 
the transiting planet; (2) the nearness to resonance of the two planets; (3) the eccentricities of the planets. The detectability 
of such planets also depends on the measurement error, the intrinsic noise due to stellar variability, and the number of transit 
timing measurements. One requirement for the case of an external perturbing planet is that observations should be made over 
a time longer than the period of the timing variations, which can be longer than the period of the perturbing planet. Ignoring 
these complications, a rough estimate of detectability can be obtained from comparing the standard deviation of the transit 
timing with the measurement error. 

It is worthwhile to provide a numerical example for the case of a hot Jupiter with a 3 day period that is perturbed by 
a lighter, exterior planet on a circular orbit in exact 2:1 resonance. The timing deviation amplitude is of order the period (3 
days) times the mass ratio (300) or about 3 minutes (equation I33I I: 

5t = 3 mj^A (3^) 

These variations accumulate over a time-scale of order the period (3 days) times the planet to star mass ratio to the power 
of 2/3, which for a transiting planet of order a Jupiter mass is about 5 months feauation l34ll : 

tcycie = 150 f ) days. (38) 
V mj J 

Such a large signal should easily be detectable from the ground. With relative photometric precision of lO"'^ from space or 



© 0000 RAS, MNRAS 000, 000-000 



Transit timing variations 13 



100 



50 - 







-50 



-100 



— 1 1 r- 



— 1 1 r- 



A 



A 

• \ • .fx 



—I 1 r- 



A * 



•./•• • 



a: 
: V 



• 

V 



_l I l_ 



_l I l_ 



—I I l_ 



50 100 

Eclipse number (P=30.1d) 



150 



200 



30 



20 



10 



— 1 1 r- 



If 



-10 



-20 



-30 



_l I l_ 



20 40 60 

Eclipse number (P=61.0cl) 



80 



100 



Figure 6. (a) Transit timing variations for Gliese 876 B and (b) Gliese 876 C. The vertical axis is in units of hours, while the horizontal 
axis is in units of the period of the transiting planet (given in parentheses for reference). 



from future ground-based experiments, less massive objects or objects further away from resonance could be detected. The 
observations could be scheduled in advance and require a modest amount of observing time with the possible payoff of being 
able to detect a terrestrial-sized planet. 



8.2 Comparison to other terrestrial planet search techniques 

To attempt a comparison with other transit timing techniques, we have estimated the mass of a planet which may bo detected 
at an amplitude of 10 times the noise for a given technique. We compare three techniques for measuring the mass of planets: 
(1) radial velocity variations of the star; (2) astrometric measurements; (3) transit timing variations (TTV). We assume that 
radial velocity measurements have a limit of 0.5 m/s, which is about the highest accuracy that has been achieved from the 
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perl' HD2094S8B 

Figure 7. Mass sensitivity of various planet detection techniques to secondary planets in HD 209458. The vertical axis is the perturbing 
planet's mass in units of Mq. The horizontal axis is the period ratio of the planets. The solid line is for the transit timing technique, the 
dashed line is astrometric, and the dotted line is the radial velocity technique. The large dots are equation 1331 . 



ground, and may be at the limit imposed by stellar variability IIButler et al J|200J). We assume that astrometric measurements 
have an accuracy of 1 /jarcsecond which is th e accuracy which is projected to be achieved by the Space Interferometry Mission 
iPord fc Tremaiii3l2003l : ISozzetti et alJl2003ll . Finally, we assume that the transit timing can be measured t o an accuracy of 
10 seconds, which is the highest accuracy of transit timing measurements of HD 209458 jBrown et al.ll200lh . 

We concentrate on HD 209458 since it is the best studied transiting planet. This system is at a distance of 46 pc and has a 
period of 3.5 days. Figure [7| shows a comparison of the 10-cr sensitivity of these three techniques. The solid curve is computed 
for mtrans = 6.7 X 1O~''M0 and nipert ~ 10~^Mq and both planets on circular orbits. When not trapped in resonance, the 
amplitude of the timing variations scales as nipert /mo, so we scale the results to the mass of the perturber to compute where 
the timing variations are 100 seconds - this determines the sensitivity. For exact resonance we plot equation 11331 . The TTV 
technique is more sensitive than the astrometric technique at semi-major axis ratios smaller than about 2. Off resonance, 
radial velocity measurements are the technique of choice for this system, while on resonance the TTV is sensitive to much 
smaller planet masses. Note that in Figure |7| the TTV and astrometric techniques have the same slope at small Ppert/ Ptrans- 
This is because the transit timing technique is measuring the reflex motion of the host star due to the inner planet, which is 
also being measured by astrometry. The solid curve is an upper limit to the minimum mass detectable in HD 209458 since a 
non-zero eccentricity will lead to larger timing variations (Figure |5J and thus a smaller detectable mass. 



8.3 Breaking the mass-radius degeneracy 

In the case that two planets are discovered to transit their host star, measurement of the transit timing variations can break 
the degener acy between mass and rad i us ne eded to derive the physical parameters of the planetary system. This has been 
discussed bv lSeaeer fc Mallen-QrnelasI (I2OO3II who use a theoretical stellar mass-radius relation to break this degeneracy. We 
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provide a simplified treatment here to illustrate the nature of the degeneracy and how it can be broken with observations of 
transit timing variations. 

Consider a planetary system with two transiting planets on circular orbits which are coplanar, exactly edge-on, and have 
measured radial velocity amplitudes. We'll assume that the star is uniform in surface brightness and that mi,m2 <^ mo. We'll 
also assume that the unperturbed periods Pi , P2 can be measured from the duration between transits. Then there are eight 
physical parameters of interest which describe the system: mo,mi,m2, Ro, Ri, R2,ai, and 02 where Ri are the radii of the 
star and planets. Without measuring the transit timing variations, there are a total often parameters which can be measured: 
Ki, K2,tTi,tT2,tgi,tg2, ^Fi, AF2, Pi, and P2, where trj labels the duration of transit from mid-ingress to mid-egress, tgj 
labels the duration of ingress or egress for planet j, Kj are the velocity amplitudes of the two planets, and AFj are the relative 
depths of the transits in units of the uneclipsed brightness of the star (for planet j) . Although there are more constraining 
parameters than model parameters, there is a degeneracy since some of the observables are redundant. All of the system 
parameters can be expressed in terms of observables and the ratio of the mass to radius of the star, mo/Ro, 



Ro ' 



1/2 1 /o 



mo Y TvtTj G J \ Ro . 

P, fGTvtTjV^^ /mo\i/2 
' 2n[ P, J \Ro) 
Ro _ ''^irj 

Rj _ TTfgj 



(39) 



wher e j — 1,2 labels each planet. From this information alone one can constrain the density of the star iSeaeer fc Mallen-Ornela^ 
boost) . For the simplified case discussed here. 



for either planet (this differs sligthly from the expression in lSeaeer fc Mallen-Ornelasl2003l . since we define the transit duration 
from mid in/egress). If, in addition, one can measure the amplitude of the transit timing variations of the outer planet, o"2, 
then this determines the mass ratio. For the case that the star's motion dominates the transit timing, 

mi _ 2^710-2 

mo"pV3p2/3- ^ > 

For other cases, the transit timing amplitude can be computed numerically. Then, from the above expression for nii/mo one 
can find the ratio of the mass to the radius of the star 

1 r>7/3 o2/3 t^2 

mo ^ 1 Pi P2' Ki 

Ro Stt^G tTiol ■ ' 

Combined with the measurement of the density, this gives the absolute mass and radius of the star. This procedure requires 
no assumptions about the mass-radius relation for the host star, and in principle could be used to measure this relation. If 
one can also measure transit timing variations for the inner planet, then an extra constraint can be obtained 

Pim2 J., , , 
ai = — — / a , 43 
mo 

where /(a) is a function derived from averaging eauation IA7I (Note that the phase of the orbits is needed for this equation, 
which can be found from the velocity amplitude curve). This provides an extra constraint on the system, and thus will be a 
check that this procedure is robust. 

Clearly we have made some drastically simplifying assumptions which are not true for any physical transit. The inclination 
of the orbits must be solve d for, which can be done f r om th e ratio of the durations of the ingress and transit and the change 
in fiux, as demonstrated bv lSeaeer fc Mallen-OrnelasI (I2OO3II. In addi t ion, li mb-darkening must be included, and can be solved 
for with high signal-to-noise data as demonstrated bv lBrown et al.l l)200lh . Finally, the orbits are not generally circular, so 
the parameters ej,zuj,Q.j,(jj, which can be derived from the velocity amplitude measurements, should be accounted for. The 
general solution is rather complicated and would best be accomplished numerically, but the degeneracy has a similar nature 
to the circular case and can in principle be broken by the transit timing variations. 
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9 EFFECTS WE HAVE IGNORED 

We now discuss several physical effects that we have ignored, which ought to be kept in mind by observers measuring transit 
timing variations. 

9.1 Light travel time 

et al carried out a search for perturbing planets in the eclipsing binary stellar system, CM Draconis, using the 

changes in the times of the eclipse due to the light travel time to measure a tentative signal consistent with a Jupiter-mass 
planet at ~ 1 AU (their technique would in principle be sensitive to a planet on an eccentric orbit as well, c.f equation I24II . 
The "Romer Effect" due to the change in light travel time caused by the reffex motion of the inner binary is much smaller in 
planetary systems than in binary stars since their masses and semi-major axes are small, having an amplitude 

where Mj is the mass of Jupiter and 02 is the semi-major axis of the perturbing planet. This effect is present in the absence 
of deviations from a Keplerian orbit because the inner binary orbits about about the center of mass. 

There can also be changes in the timing of the transit as the distance of the transiting planet from the star varies. In 
this case, the time of transit is delayed by the light travel time between the different locations where the planet intercepts 
the beam of light from the star. The amplitude of these variations is smaller than the o we have calculated by a factor of 
~ Vtrans/c, whcrc vtrans IS the velocity of the transiting planet. So, only very precise measurements will require taking into 
account light travel time effects, which should be borne in mind in future experiments (of course the light-travel time due to 
the motion of the observer in our solar system must be taken into account with current experiments). 



9.2 Inclination 

We have assumed that the planets are strictly coplanar and exactly edge-on. The first assumption is based on the fact that 
the solar system is nearly coplanar and the theoretical pre judice that plan ets forming out of disks should be nearly coplanar. 
Small non-coplanar effects will change our results slightly iMiralda-Escudeii2002') . while large inclination effects would require 
a reworking of the theory. Since some extrasolar planetary systems have been found with rather large eccentricities, it is 
entirely possible that non-coplanar systems will be found as well, a possibility we leave for future studies. 

The assumption that the systems are edge-on is based on the fact that a transit can occur only for systems that are 
nearly edge-on. For small inclinations our formulae will only change slightly, but may result in interesting effects such as a 
change in the duration of a transit, or even the disappearance of transits due to the motion of the star about the barycenter 
of the system. On a much longer time-scale (centuries), the precession of an eccentric orbit might cause th e disappearance 
of tra nsits since the projected shape of the orbit on the sky can change. This possibility was mentioned bv lLauehlin et all 
I^2nn4^ for GJ 876. 



9.3 Other sources of timing "noise" 

Aside from the long term effects that have been ignored there are several sources of timing noise that must be included in 
the analysis of observations of transiting systems. These sources of noise could come from the planet or the host star. If 
the planet has a moon or is a binary planet then there is some wobble in its position causing a change in both the timing 
and duration of a transit JSartoretti fc Schneideilfl999l : Isrown et al.ii200L) . A moon or ring s ystem may transit before the 
planet causing a shallower tra nsit to appear earlier or later than it would without the moon dSrown et aljEoOll : ISchneideil 
1200,'j: iBarnes fc Fortnell2004h . A large scale asymmetry of the planet's shape with respect to its center of mass might cause 
a slightly earlier or later start to the ingress or end of egress. 

Stellar variability could also make a significant contribution to the noise. Variations in the brightness of the star might 
affect the accuracy of the measurement of the start of ingress and the end of egress, which are the times that are critical to 
timing of a transit. Stellar oscillations can cause variations in the surface of the Sun of ~ 100 km in regions of size 10'^~* km, 
which corresponds to a one second variation for a planet moving at 100 km/s. 



9.4 Coverage Gaps 

Radial velocities and prior transit lightcurves predict the epoch of future transits and to schedule photometric monitoring for 
the system of interest. Observational limitations (e.g. bad weather, equipment failure, scheduling requirements) will lead to 
transits being missed, which in turn will cause inaccuracies in a. Since the signal is periodic, 5t{t) may be straightforward 
to extract with a few missed transits; however, if the outer planet is highly eccentric, then most of the change in transit 
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timitiK may occ ur for a few transits (e.g. HD37124; a similar selection effect occurs in radial yelocity searches as discussed by 
ICumminell200l . In principle this effect will ayerage out oyer long obseryational interyals; howeyer as in this context "long" 
may mean seyeral decades or more, it will be important to eyaluate the effect of coyerage gaps on detections oyer a time-scale 
of months-years. We will return to this in detail in future work. 

We note in passing that the adyent of the new astrometric all-sky suryeys such as Gaia ijPerryman et all2nnitl will proyide 
photometry for ~ lO'' stars > 15 mag, with fewer coyerage gaps than ground-based obseryations; we thus expect the detection 
method by transits alone ( section |8.H to really come into its own oyer the next two decades. Assuming ~ 0.4 detections (three 
transits) per 10* stars (c.f. iBrowrjboO Sl . we may expect transit lightcuryes of perhaps ~ 1000 exoplanetary systems oyer the 
mission lifetime, greatly aiding the determinatio n of a for these syste ms. The NASA Kepler mission will also proyide uniform 
monitoring of about 10'^"'^ transiting gas giants jBorucki et aljbooall . a nd if flown, the Microle nsing Planet Finder (formerly 
GEST) will discoyer 10*~^ transiting gas giants with uniform coyerage jBennett fc et al.ll2003l) . 



10 CONCLUSIONS 

For an exoplanetary system where one or more planets transit the host star the timing and duration of the transits can be 
used to deriye seyeral physical characteristics of the system. This technique breaks the degeneracy between the mass and 
radius of the objects in the system. The inclination, absolute mass, and absolute radii of the star and planets can be found; 
in principle this could be used to measure the mass-radius relation for stars that are not in eclipsing binaries. 

In addition, TTV can be used to infer the existence of preyiously undetected planets. We haye found that for yariations 
which occur oyer seyeral orbital periods the strongest signals occur when the perturbing planet is either in a mean-motion 
resonance with the transiting planet or if the transiting planet has a long period (which, unfortunately, makes a transit less 
probable). The resonant case is more interesting since the probability of a planet transiting decreases significantly as the 
semi-major axis becomes large. Using the TTV scheme it is possible to detect earth-mass planets using current obseryational 
technologies for both ground based and space based obseryatories. Data for seyeral transits of HD 209458 could be gathered 
and studied oyer a relatiyely short time due to its small period. Once the existence of a second planet is established one can 
predict the times at which it would likely transit the host star. Follow-up obseryations with HST or high-precision ground-based 
telescopes at those times would increase the likelihood of detecting a transit of the second planet. 

If the second planet is terrestrial in nature, this transit timing method may be the only way currently to determine the 
mass of such planets in other star systems. Astrometry is another possible technique but it may take a decade of technological 
deyelopment before the necessary sensitiyity is achieyed. In addition, complementary techniques are necessary to probe 
different parts of parameter spac e and to proyide extra confidence that the detected planets are real, giyen the likely low 
signal-to-noise iCould et al.ll2'o04h . For the near future the TTV technique may be the most promising method of detecting 
earth-mass planets around main sequence stars besides the Sun. 

We exhort obseryers to (1) discoyer longer period transiting planets JSeagroyes et al.l l2003*l since the signal increases 
with transiting planet's period; (2) increase the signal-to-noise of ground based differential photometry ijllowell et alJl200^ 
for more precise m easurement of the transit times; and (3) examine their transit data for the presence of perturbing planets 
ijBrown et al.ll200ll) . 

The treatment of this problem has ignored many effects which we plan to take into account in a future paper in which 
we will simulate realistic lightcuryes including noise and to fit the simulated data to deriye the parameters of the perturbing 
planet, exploring degeneracies in the period ratio. We will also deriye the probability of detecting such systems taking into 
account yarious assumptions about the formation, eyolution, and stability of extrasolar planets. 
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APPENDIX A: PERTURBATION THEORY TREATMENT OF NEARLY CIRCULAR ORBITS 

In this appendix we consider the case of two planets whose orbits are nearly circular and whose timing yariations are dominated 
by cha nges in eccentricity (see • The timing yariations can be computed from a Hamiltonian as described in iMalhotral 
Jl993ah . We keep terms which are first order in the eccentricity because mutual perturbations between the planets induces 
an eccentricity of order 1112/ mo- To first order in the eccentricities, the Taylor-expanded Hamiltonian is^ 

H^nt = -^1^ [2P-2acosV'-eiC+cos(Ai -tui) 

^ We have corrected equation (26) in lMalhotral ll993a^ which should have a — ocosi/i in the second line. 



© 0000 RAS, MNRAS 000, 000-000 



18 Agol et al. 



+ e2do cos (A2 — 072) 

ei {c~ + a5j,x) cos ((j + l)Ai - j\2 - roi) 

ei (cj+ - 3q(5j,i) cos ((j - l)Ai - jA2 + roi) 

+ 62 (4 - 4a<S,- 1) cos (jAi - (j + 1)A2 + 1172) 



(Al) 



where Si^j is the Kronecker deha, = dab^^j,^ ^ ^-^^i/lj' ^i/2('-*)' ~ 11/12, feJ^^jC'^) ^^e Laplace coefficient, 

cos 



/(, Vl — 2a cos 5 + 
where ■;/; = Ai — A2 , 

2\-l/2 



and 



2q cos ijj + a 



d 



This equation includes no secular terms since these are higher order in the eccentricity. Note that since we have only included 
the first order terms in the eccentricity, the resonant arguments which appear have ratios j + l:j and j:j + 1 for the mean 
longitudes. 

The perturbed semi-major axis is given in lMalhotraUl993bll . and we compute the perturbed eccentricity and longitude of 
periastron using h\ — ei sin^zri and ki = ei cos^zri. Keeping all the resonance terms that exist to first order in the eccentricities 
gives the equations of motion for h\,k\, 



(A2) 



(A3) 



(A4) 



hi 



1 ^2 r + 

— n\a Cn cosAi 

9 L 



+ 
+ 



fcl 



mo 

00 

^ { (c+ - 3a<5,,i) cos [{j - l)Ai - jA2] 
i=i 

(c~ + Q(5j,i) COS [(j + l)Ai - jA2]}] , 

1 m2 r + . 

-niQ c,, smAi 

2 mo '- 



E 



{ (cj^ ~ SofJj,!) sin [(j - l)Ai - jA2] 



(A5) 



(c^. +a<5i,i)sin[(j + l)Ai - jA2]}] 



To find the change in the transit timing we use the orbital elements to compute the variation in 9\. To first order in ei 
5n\ 



— 2- ~ 1 H + 2fci cos [nu)t + Aio] + 2h\ sin [nioi + A 

mo riio 



(A6) 



Since we begin with zero eccentricity, we ignore perturbations to A in the sin and cos terms in this equation. As in equation 
1221 1 where 56\ = 9i — nio , we integrate this equation to find 



5ti 



3 mo 
— -a 



(q(^) 



24iK{k) 



ct sin x}} 



mo 



+ 



sin (Al — Aio) 

"1 

Eoo gj -Sq^j,! / sin [jj,] _ sin [nit+j^g] 
j — 1 (j — l)ni-jn2 \^j(ni— n 

y^oo "7 +"'^3.1 / sin 
= l (i + l)"l-i"-2 



2) "1 

[i'/'] _ sin [nit-ji/»o] 
-"2) 



(A7) 



where ni and 712 are taken at their initial values, A io = Ai(t = 0), A2 = A2(t = 0), k = 2^/a/{l + a), K{k) is the complete 
elliptic integral, Q{'4>) is defined in the appendix of iMalhotral il993lJ) . and we have dropped any terms which vary linearly 
with time. 

A similar calculation can be carried out for perturbations by a planet interior to the transiting planet. 



5t2 = 



-3i^ (TTT^ (O(V') " - «sin(^)) 



722- 



l[A2 — A2 
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Figure Al. Deviations from uniform times between transits for a transiting outer planet with semi-major axis ratio 1.8. Crosses represent 
the analytic result, equation IA8i . while the diamonds are from a numerical integration of the equations of motion. The horizontal axis 
is the number of the transit, while the vertical axis shows the timing differences in units of Pi, the period of the inner perturbing planet, 
multiplied by mo/mpert- 



Z^j = l (ini-(j + l)n2) yj(ni-n 
j = l (jni-(j-l)n2) yj(ni-n2 



Bin[n2t+jiio] 
"2 

3in["2t-jV)o] 
"2 



(A8) 



A comparison of these equations to numerical calculations is shown in Figure IXTl The planets are on initially circular orbits 
with a semi-major axis ratio of 1.8, or period ratio of 2.4. The masses are equal and the planets start aligned along the line 
of sight at the first transit. 

In the case that the semi-major axis dominates the timing variations, one can take perturbed value of the eccentricity 
{h and k) and compute the change in semi-major axis due to these eccentricities. The result is a double-series over resonant 
terms, so we do not reproduce it here. 



REFERENCES 

Alonso R., Brown T. M., Torres G., Latham D. W., Sozzetti A., Mandushev G., Belmonte J. A., Charbonneau D., Deeg 

H. J., Dunham E. W., O'Donovan F. T., Stefanik R. R, 2004, astro-ph/0408421 
Barnes J. W., Fortney J. J., 2004, astro-ph/0409506 

Bennett D. P., et al. 2003, in Future EUV/UV and Visible Space Astrophysics Missions and Instrumentation. Edited by J. 

Chris Blades, Oswald H. W. Siegmund. Proceedings of the SPIE, Volume 4854, pp. 141-157 (2003). The Galactic Exoplanet 

Survey Telescope (GEST). pp 141-157 
Borkovits T., Erdi B., Forgacs-Dajka E., Kovacs T., 2003, A&A, 398, 1091 



© 0000 RAS, MNRAS 000, 000-000 



20 Agol et al. 



Borucki W. J., Koch D. G., Lissauer J. J., Basri G. B., Caldwell J. F., Cochran W. D., Dunham E. W., Geary J. C, Latham 
D. W., Gilliland R. L., Caldwell D. A., Jenkins J. M., Kondo Y., 2003, in Future EUV/UV and Visible Space Astrophysics 

Missions and Instrumentation. Edited by J. Chris Blades, Oswald H. W. Siogmund. Proceedings of the SPIE, Volume 4854, 
pp. 129-140 (2003). The Kepler mission: a wide-field-of-view photometer designed to determine the frequency of Earth-size 
planets around solar-like stars, pp 129-140 

Bouchy F., Pont F., Santos N. C, Mclo C, Mayor M., Queloz D., Udry S., 2004, astro-ph/0404264 

Brown T., Charbonneau D., Gilliland R. L., Noyes R. W., A. B., 2001, ApJ, 552, 699 

Brown T. M., 2003, ApJL, 593, L125 

Butler P., Vogt S. S., Marcy G. W., Fischer D. A., Wright J. T., Henry G. W., Laughlin G., Lissauer J., 2004, astro-ph/0408587 
Butler R. P., Marcy G. W., Fischer D. A., Brown T. M., Contos A. R., Korzennik S. G., Nisenson P., Noyes R. W., 1999, 
ApJ, 526, 916 

Charbonneau D., 2003, astro-ph/0312252 

Charbonneau D., Brown T., Latham D., Mayor M., 2000, ApJ, 529, L45 
Gumming A., 2004, astro-ph/0408470 

Deeg H. J., Doyle L. R., Kozhevnikov V. P., Blue J. E., Martin E. L., Schneider J., 2000, A&A, 358, L5 

Dobrovolskis A. R., Borucki W. J., 1996, American Astronomical Society Division of Planetary Science Meeting, 28, 1208 

Fischer D. A., Marcy G. W., Butler R. P., Vogt S. S., Frink S., Apps K., 2001, ApJ, 551, 1107 

Ford E. B., Tremaine S., 2003, PASP, 115 

Gould A., Gaudi B. S., Han C, 2004, ArXiv Astrophysics e-prints 

Horne K., 2003, in ASP Conf. Ser. 294: Scientific Frontiers in Research on Extrasolar Planets Status aand Prospects of 

Planetary Transit Searches: Hot Jupiters Galore 
Howell S. B., Everett M. E., Tonry J. L., Pickles A., Dain C, 2003, PASP, 115, 1340 
Konacki M., Torres G., Jha S., Sasselov D. D., 2003, Nature, 421, 507 

Konacki M., Torres G., Sasselov D. D., Pietrzynski G., Udalski A., Jha S., Ruiz M. T., Gieren W., Minniti D., 2004, 
astro-ph/0404541 

Laughlin G., Butler R. P., Fischer D. A., Marcy G. W., Vogt S. S., Wolf A. S., 2004, astro-ph/0407441 

Laughlin G., Chambers J. E., 2001, ApJ, 551, L109 

Malhotra R., 1993a, in Planets Around Pulsars, ASP Conference Proceedings, J. A. Phillips, J. E. Thorsett, and S. R. 

Kulkaxni, eds., p. 89 
Malhotra R., 1993b, ApJ, 407, 266 
Mandel K., Agol E., 2002, ApJ, 529, L45 

Marcy G. W., Butler R. P., Fischer D., Vogt S. S., Lissauer J. J., Rivera E. J., 2001, ApJ, 556, 296 

Marcy G. W., Butler R. P., Fischer D. A., Laughhn G., Vogt S. S., Henry G. W., Pourbaix D., 2002, ApJ, 581, 1375 
McArthur B. E., Endl M., Cochran W. D., Benedict G. F., Fischer D. A., Marcy G. W., Butler R. P., Naef D., Mayor M., 

Queloz D., Udry S., Harrison T. E., 2004, astro-ph/0408585 
Miralda-Escudo J., 2002, ApJ, 564, 1019 

Murray C. D., Dermott S. F., 1999, Solar System Dynamics. (Cambridge: Cambridge Univ. Press) 

Nauenberg M., 2002, ApJ, 568, 369 

Pan M., Sari R., 2004, A J, 128, 1418 

Ferryman M. A. C, 2000, Rep. Prog. Phys., 63, 1209 

Ferryman M. A. C, de Boer K. S., Gilmore G., H0g E., Lattanzi M. G., Lindegren L., Luri X., Mignard F., Pace O., de 

Zeeuw P. T., 2001, A&A, 369, 339 
Pont F., Bouchy F., Queloz D., Santos N., Melo C, Mayor M., Udry S., 2004, astro- ph/0408499 
Rivera E. J., Lissauer J. J., 2001, ApJ, 558, 392 

Santos N. C, Bouchy F., Mayor M., Pepe F., Queloz D., Udry S., Lovis C, Bazot M., Benz W., Bertaux J. ., Curto G. L., 

Delfosse X., Mordasini C, Naef D., Sivan J. ., Vauclair S., 2004, astro-ph/0408471 
Sartoretti P., Schneider J., 1999, A&AS, 134, 553 

Schneider J., 2003, SF2A-2003: Semaine de I'Astrophysique Frangaise, meeting held in Bordeaux, France, June 16-20, 2003. 

Eds.: F. Combes, D. Barret and T. Contini. EdP-Sciences, Conference Series, p. 75 
Seager S., Mallcu-Ornelas G., 2003, ApJ, 585 

Seagroves S., Harker J., Laughlin G., Lacy J., Castellano T., 2003, PASP, 115, 1355 
Sozzetti A., Casertano S., Brown R. A., Lattanzi M. G., 2003, PASP, 115 
Vogt S. S., Marcy G. W., Butler R. P., Apps K., 2000, ApJ, 536, 902 
Wisdom J., 1980, AJ, 85, 1122 



© 0000 RAS, MNRAS 000, 000-000 



